<!DOCTYPE HTML>
<html>
<head>
<meta charset="UTF-8">
<title>FDLA and FMMC solutions for a 50-node, 200-edge graph</title>
<link rel="canonical" href="http://cvxr.com/cvx/examples/graph_laplacian/html/larger_example.html">
<link rel="stylesheet" href="../../examples.css" type="text/css">
</head>
<body>
<div id="header">
<h1>FDLA and FMMC solutions for a 50-node, 200-edge graph</h1>
Jump to:&nbsp;&nbsp;&nbsp;&nbsp;
<a href="#source">Source code</a>&nbsp;&nbsp;&nbsp;&nbsp;
<a href="#output">Text output</a>
&nbsp;&nbsp;&nbsp;&nbsp;
<a href="#plots">Plots</a>
&nbsp;&nbsp;&nbsp;&nbsp;<a href="../../index.html">Library index</a>
</div>
<div id="content">
<a id="source"></a>
<pre class="codeinput">
<span class="comment">% S. Boyd, et. al., "Convex Optimization of Graph Laplacian Eigenvalues"</span>
<span class="comment">% ICM'06 talk examples (www.stanford.edu/~boyd/cvx_opt_graph_lapl_eigs.html)</span>
<span class="comment">% Written for CVX by Almir Mutapcic 08/29/06</span>
<span class="comment">% (figures are generated)</span>
<span class="comment">%</span>
<span class="comment">% In this example we consider a graph described by the incidence matrix A.</span>
<span class="comment">% Each edge has a weight W_i, and we optimize various functions of the</span>
<span class="comment">% edge weights as described in the referenced paper; in particular,</span>
<span class="comment">%</span>
<span class="comment">% - the fastest distributed linear averaging (FDLA) problem (fdla.m)</span>
<span class="comment">% - the fastest mixing Markov chain (FMMC) problem (fmmc.m)</span>
<span class="comment">%</span>
<span class="comment">% Then we compare these solutions to the heuristics listed below:</span>
<span class="comment">%</span>
<span class="comment">% - maximum-degree heuristic (max_deg.m)</span>
<span class="comment">% - constant weights that yield fastest averaging (best_const.m)</span>
<span class="comment">% - Metropolis-Hastings heuristic (mh.m)</span>

<span class="comment">% randomly generate a graph with 50 nodes and 200 edges</span>
<span class="comment">% and make it pretty for plotting</span>
n = 50; threshold = 0.2529;
rand(<span class="string">'state'</span>,209);
xy = rand(n,2);

angle = 10*pi/180;
Rotate = [ cos(angle) sin(angle); -sin(angle) cos(angle) ];
xy = (Rotate*xy')';

Dist = zeros(n,n);
<span class="keyword">for</span> i=1:(n-1);
  <span class="keyword">for</span> j=i+1:n;
    Dist(i,j) = norm( xy(i,:) - xy(j,:) );
  <span class="keyword">end</span>;
<span class="keyword">end</span>;
Dist = Dist + Dist';
Ad = Dist &lt; threshold;
Ad = Ad - eye(n);
m = sum(sum(Ad))/2;

<span class="comment">% find the incidence matrix</span>
A = zeros(n,m);
l = 0;
<span class="keyword">for</span> i=1:(n-1);
  <span class="keyword">for</span> j=i+1:n;
    <span class="keyword">if</span> Ad(i,j)&gt;0.5
      l = l + 1;
      A(i,l) =  1;
      A(j,l) = -1;
    <span class="keyword">end</span>;
  <span class="keyword">end</span>;
<span class="keyword">end</span>;
A = sparse(A);

<span class="comment">% Compute edge weights: some optimal, some based on heuristics</span>
[n,m] = size(A);

[ w_fdla, rho_fdla ] = fdla(A);
[ w_fmmc, rho_fmmc ] = fmmc(A);
[ w_md,   rho_md   ] = max_deg(A);
[ w_bc,   rho_bc   ] = best_const(A);
[ w_mh,   rho_mh   ] = mh(A);

tau_fdla = 1/log(1/rho_fdla);
tau_fmmc = 1/log(1/rho_fmmc);
tau_md   = 1/log(1/rho_md);
tau_bc   = 1/log(1/rho_bc);
tau_mh   = 1/log(1/rho_mh);

eig_opt  = sort(eig(eye(n) - A * diag(w_fdla) * A'));
eig_fmmc = sort(eig(eye(n) - A * diag(w_fmmc) * A'));
eig_mh   = sort(eig(eye(n) - A * diag(w_mh)   * A'));
eig_md   = sort(eig(eye(n) - A * diag(w_md)   * A'));
eig_bc   = sort(eig(eye(n) - A * diag(w_bc)   * A'));

fprintf(1,<span class="string">'\nResults:\n'</span>);
fprintf(1,<span class="string">'FDLA weights:\t\t rho = %5.4f \t tau = %5.4f\n'</span>,rho_fdla,tau_fdla);
fprintf(1,<span class="string">'FMMC weights:\t\t rho = %5.4f \t tau = %5.4f\n'</span>,rho_fmmc,tau_fmmc);
fprintf(1,<span class="string">'M-H weights:\t\t rho = %5.4f \t tau = %5.4f\n'</span>,rho_mh,tau_mh);
fprintf(1,<span class="string">'MAX_DEG weights:\t rho = %5.4f \t tau = %5.4f\n'</span>,rho_md,tau_md);
fprintf(1,<span class="string">'BEST_CONST weights:\t rho = %5.4f \t tau = %5.4f\n'</span>,rho_bc,tau_bc);

<span class="comment">% plot results</span>
figure(1), clf
gplot(Ad,xy);
hold <span class="string">on</span>;
plot(xy(:,1), xy(:,2), <span class="string">'ko'</span>,<span class="string">'LineWidth'</span>,4, <span class="string">'MarkerSize'</span>,4);
axis([0.05 1.1 -0.1 0.95]);
title(<span class="string">'Graph'</span>)
hold <span class="string">off</span>;

figure(2), clf
v_fdla = [w_fdla; diag(eye(n) - A*diag(w_fdla)*A')];
[ifdla, jfdla, neg_fdla] = find( v_fdla.*(v_fdla &lt; -0.001 ) );
v_fdla(ifdla) = [];
wbins = [-0.6:0.012:0.6];
hist(neg_fdla,wbins); hold <span class="string">on</span>,
h = findobj(gca,<span class="string">'Type'</span>,<span class="string">'patch'</span>);
set(h,<span class="string">'FaceColor'</span>,<span class="string">'r'</span>)
hist(v_fdla,wbins); hold <span class="string">off</span>,
axis([-0.6 0.6 0 12]);
xlabel(<span class="string">'optimal FDLA weights'</span>);
ylabel(<span class="string">'histogram'</span>);

figure(3), clf
xbins = (-1:0.015:1)';
ymax  = 6;
subplot(3,1,1)
hist(eig_md, xbins); hold <span class="string">on</span>;
max_md = max(abs(eig_md(1:n-1)));
plot([-max_md -max_md],[0 ymax], <span class="string">'b--'</span>);
plot([ max_md  max_md],[0 ymax], <span class="string">'b--'</span>);
axis([-1 1 0 ymax]);
text(0,5,<span class="string">'MAX DEG'</span>);
title(<span class="string">'Eigenvalue distributions'</span>)
subplot(3,1,2)
hist(eig_bc, xbins); hold <span class="string">on</span>;
max_opt = max(abs(eig_bc(1:n-1)));
plot([-max_opt -max_opt],[0 ymax], <span class="string">'b--'</span>);
plot([ max_opt  max_opt],[0 ymax], <span class="string">'b--'</span>);
axis([-1 1 0 ymax]);
text(0,5,<span class="string">'BEST CONST'</span>);
subplot(3,1,3)
hist(eig_opt, xbins); hold <span class="string">on</span>;
max_opt = max(abs(eig_opt(1:n-1)));
plot([-max_opt -max_opt],[0 ymax], <span class="string">'b--'</span>);
plot([ max_opt  max_opt],[0 ymax], <span class="string">'b--'</span>);
axis([-1 1 0 ymax]);
text(0,5,<span class="string">'FDLA'</span>);

figure(4), clf
xbins = (-1:0.015:1)';
ymax  = 6;
subplot(3,1,1)
hist(eig_md, xbins); hold <span class="string">on</span>;
max_md = max(abs(eig_md(1:n-1)));
plot([-max_md -max_md],[0 ymax], <span class="string">'b--'</span>);
plot([ max_md  max_md],[0 ymax], <span class="string">'b--'</span>);
axis([-1 1 0 ymax]);
text(0,5,<span class="string">'MAX DEG'</span>);
title(<span class="string">'Eigenvalue distributions'</span>)
subplot(3,1,2)
hist(eig_mh, xbins); hold <span class="string">on</span>;
max_opt = max(abs(eig_mh(1:n-1)));
plot([-max_opt -max_opt],[0 ymax], <span class="string">'b--'</span>);
plot([ max_opt  max_opt],[0 ymax], <span class="string">'b--'</span>);
axis([-1 1 0 ymax]);
text(0,5,<span class="string">'MH'</span>);
subplot(3,1,3)
hist(eig_fmmc, xbins); hold <span class="string">on</span>;
max_opt = max(abs(eig_fmmc(1:n-1)));
plot([-max_opt -max_opt],[0 ymax], <span class="string">'b--'</span>);
plot([ max_opt  max_opt],[0 ymax], <span class="string">'b--'</span>);
axis([-1 1 0 ymax]);
text(0,5,<span class="string">'FMMC'</span>);

figure(5), clf
v_fmmc = [w_fmmc; diag(eye(n) - A*diag(w_fmmc)*A')];
[ifmmc, jfmmc, nonzero_fmmc] = find( v_fmmc.*(v_fmmc &gt; 0.001 ) );
hist(nonzero_fmmc,80);
axis([0 1 0 10]);
xlabel(<span class="string">'optimal positive FMMC weights'</span>);
ylabel(<span class="string">'histogram'</span>);

figure(6), clf
An = abs(A*diag(w_fmmc)*A');
An = (An - diag(diag(An))) &gt; 0.0001;
gplot(An,xy,<span class="string">'b-'</span>); hold <span class="string">on</span>;
h = findobj(gca,<span class="string">'Type'</span>,<span class="string">'line'</span>);
set(h,<span class="string">'LineWidth'</span>,2.5)
gplot(Ad,xy,<span class="string">'b:'</span>);
plot(xy(:,1), xy(:,2), <span class="string">'ko'</span>,<span class="string">'LineWidth'</span>,4, <span class="string">'MarkerSize'</span>,4);
axis([0.05 1.1 -0.1 0.95]);
title(<span class="string">'Subgraph with positive transition prob.'</span>)
hold <span class="string">off</span>;
</pre>
<a id="output"></a>
<pre class="codeoutput">
 
Calling sedumi: 2551 variables, 2350 equality constraints
------------------------------------------------------------
SeDuMi 1.21 by AdvOL, 2005-2008 and Jos F. Sturm, 1998-2003.
Alg = 2: xz-corrector, Adaptive Step-Differentiation, theta = 0.250, beta = 0.500
Split 1 free variables
eqs m = 2350, order n = 103, dim = 5003, blocks = 3
nnz(A) = 3203 + 0, nnz(ADA) = 3318750, nnz(L) = 1660550
 it :     b*y       gap    delta  rate   t/tP*  t/tD*   feas cg cg  prec
  0 :            2.96E-01 0.000
  1 :   3.60E+00 2.38E-02 0.000 0.0804 0.9900 0.9900  -0.22  1  1  8.8E-01
  2 :   1.12E+00 9.42E-03 0.000 0.3954 0.9000 0.9000   3.06  1  1  1.5E-01
  3 :   9.43E-01 3.44E-03 0.000 0.3655 0.9000 0.9000   1.86  1  1  4.1E-02
  4 :   9.27E-01 1.01E-03 0.000 0.2938 0.9000 0.9000   1.09  1  1  1.2E-02
  5 :   9.08E-01 3.29E-04 0.000 0.3250 0.9000 0.9000   1.06  1  1  3.9E-03
  6 :   9.04E-01 1.05E-04 0.000 0.3193 0.9000 0.9000   1.02  1  1  1.3E-03
  7 :   9.03E-01 2.58E-05 0.000 0.2460 0.9022 0.9000   1.00  1  1  3.4E-04
  8 :   9.02E-01 4.63E-06 0.000 0.1796 0.9074 0.9000   1.00  1  1  9.1E-05
  9 :   9.02E-01 5.57E-07 0.000 0.1201 0.9187 0.9000   1.00  1  1  2.3E-05
 10 :   9.02E-01 1.12E-07 0.000 0.2011 0.9160 0.9000   1.00  1  1  5.4E-06
 11 :   9.02E-01 2.56E-08 0.000 0.2284 0.9055 0.9000   1.00  1  1  1.3E-06
 12 :   9.02E-01 5.43E-09 0.000 0.2123 0.9054 0.9000   1.00  1  1  2.7E-07
 13 :   9.02E-01 1.49E-09 0.000 0.2756 0.9063 0.9000   1.00  1  1  7.6E-08
 14 :   9.02E-01 3.41E-10 0.000 0.2279 0.9027 0.9000   1.00  1  1  1.7E-08
 15 :   9.02E-01 6.99E-11 0.000 0.2051 0.9000 0.9031   1.00  1  2  3.6E-09

iter seconds digits       c*x               b*y
 15     20.0   Inf  9.0207867670e-01  9.0207869224e-01
|Ax-b| =   3.2e-09, [Ay-c]_+ =   1.5E-09, |x|=  1.1e+01, |y|=  1.2e+00

Detailed timing (sec)
   Pre          IPM          Post
3.300E-01    2.003E+01    1.000E-02    
Max-norms: ||b||=9.400000e-01, ||c|| = 1,
Cholesky |add|=0, |skip| = 0, ||L.L|| = 54.6293.
------------------------------------------------------------
Status: Solved
Optimal value (cvx_optval): +0.902079
 
Calling sedumi: 2801 variables, 2600 equality constraints
------------------------------------------------------------
SeDuMi 1.21 by AdvOL, 2005-2008 and Jos F. Sturm, 1998-2003.
Alg = 2: xz-corrector, Adaptive Step-Differentiation, theta = 0.250, beta = 0.500
Split 1 free variables
eqs m = 2600, order n = 353, dim = 5253, blocks = 3
nnz(A) = 3702 + 0, nnz(ADA) = 3918898, nnz(L) = 1961925
 it :     b*y       gap    delta  rate   t/tP*  t/tD*   feas cg cg  prec
  0 :            8.59E-02 0.000
  1 :   1.13E+00 5.63E-02 0.000 0.6559 0.9000 0.9000   5.60  1  1  1.7E+00
  2 :   1.06E+00 2.85E-02 0.000 0.5056 0.9000 0.9000   1.54  1  1  8.3E-01
  3 :   9.31E-01 8.13E-03 0.000 0.2855 0.9000 0.9000   1.80  1  1  1.6E-01
  4 :   9.53E-01 1.71E-03 0.000 0.2103 0.9000 0.9000   1.48  1  1  2.7E-02
  5 :   9.39E-01 4.28E-04 0.000 0.2503 0.9000 0.9000   1.12  1  1  6.5E-03
  6 :   9.29E-01 2.20E-04 0.000 0.5135 0.9000 0.9000   1.06  1  1  3.3E-03
  7 :   9.24E-01 1.24E-04 0.000 0.5648 0.9000 0.9000   1.04  1  1  1.9E-03
  8 :   9.24E-01 3.10E-05 0.383 0.2500 0.9000 0.0000   1.03  1  1  1.4E-03
  9 :   9.22E-01 7.28E-06 0.012 0.2347 0.9480 0.9000   1.02  1  1  4.9E-04
 10 :   9.18E-01 2.65E-06 0.000 0.3643 0.9000 0.9071   1.02  1  1  1.8E-04
 11 :   9.17E-01 1.39E-06 0.000 0.5231 0.9255 0.9000   1.01  1  1  9.1E-05
 12 :   9.16E-01 8.29E-07 0.000 0.5969 0.9468 0.9000   1.01  1  1  5.3E-05
 13 :   9.16E-01 4.24E-07 0.225 0.5120 0.0000 0.9000   1.00  1  1  2.9E-05
 14 :   9.15E-01 1.45E-07 0.000 0.3429 0.9000 0.9098   1.00  1  1  1.0E-05
 15 :   9.15E-01 5.76E-08 0.000 0.3957 0.8138 0.9000   1.00  1  1  4.2E-06
 16 :   9.15E-01 1.96E-08 0.000 0.3407 0.9000 0.8329   1.00  1  1  1.4E-06
 17 :   9.15E-01 6.23E-09 0.000 0.3176 0.8328 0.9000   1.00  1  1  4.5E-07
 18 :   9.15E-01 1.51E-09 0.000 0.2427 0.9000 0.9000   1.00  1  1  1.1E-07
 19 :   9.15E-01 1.17E-10 0.000 0.0774 0.9900 0.9900   1.00  1  2  8.5E-09

iter seconds digits       c*x               b*y
 19     35.0   Inf  9.1515175405e-01  9.1515175583e-01
|Ax-b| =   6.9e-09, [Ay-c]_+ =   4.9E-09, |x|=  1.1e+01, |y|=  1.2e+00

Detailed timing (sec)
   Pre          IPM          Post
4.200E-01    3.498E+01    2.000E-02    
Max-norms: ||b||=9.600000e-01, ||c|| = 1,
Cholesky |add|=0, |skip| = 0, ||L.L|| = 81.924.
------------------------------------------------------------
Status: Solved
Optimal value (cvx_optval): +0.915152

Results:
FDLA weights:		 rho = 0.9021 	 tau = 9.7037
FMMC weights:		 rho = 0.9152 	 tau = 11.2784
M-H weights:		 rho = 0.9489 	 tau = 19.0839
MAX_DEG weights:	 rho = 0.9706 	 tau = 33.5236
BEST_CONST weights:	 rho = 0.9470 	 tau = 18.3549
</pre>
<a id="plots"></a>
<div id="plotoutput">
<img src="larger_example__01.png" alt=""> <img src="larger_example__02.png" alt=""> <img src="larger_example__03.png" alt=""> <img src="larger_example__04.png" alt=""> <img src="larger_example__05.png" alt=""> <img src="larger_example__06.png" alt=""> 
</div>
</div>
</body>
</html>